home *** CD-ROM | disk | FTP | other *** search
/ EnigmA Amiga Run 1999 March / EnigmA AMIGA RUN 35 (1999)(G.R. Edizioni)(IT)[!][issue 1999-03].iso / earcd / grafica / amhelios / ray_cast.cpp < prev    next >
C/C++ Source or Header  |  1999-01-01  |  9KB  |  344 lines
  1. ////////////////////////////////////////////////////////////
  2. //
  3. //  RAY_CAST.CPP - Ray Cast Form Factor Class
  4. //
  5. //  Version:    1.03A
  6. //
  7. //  History:    94/08/23 - Version 1.00A release.
  8. //              94/12/16 - Version 1.01A release.
  9. //              95/02/05 - Version 1.02A release.
  10. //              95/07/21 - Version 1.02B release.
  11. //              96/02/14 - Version 1.02C release.
  12. //              96/04/01 - Version 1.03A release.
  13. //
  14. //  Compilers:  Microsoft Visual C/C++ Professional V1.5
  15. //              Borland C++ Version 4.5
  16. //
  17. //  Author:     Ian Ashdown, P.Eng.
  18. //              byHeart Software Limited
  19. //              620 Ballantree Road
  20. //              West Vancouver, B.C.
  21. //              Canada V7S 1W3
  22. //              Tel. (604) 922-6148
  23. //              Fax. (604) 987-7621
  24. //
  25. //  Copyright 1994-1996 byHeart Software Limited
  26. //
  27. //  The following source code has been derived from:
  28. //
  29. //    Ashdown, I. 1994. Radiosity: A Programmer's
  30. //    Perspective. New York, NY: John Wiley & Sons.
  31. //
  32. //  It may be freely copied, redistributed, and/or modified
  33. //  for personal use ONLY, as long as the copyright notice
  34. //  is included with all source code files.
  35. //
  36. ////////////////////////////////////////////////////////////
  37.  
  38. #include "ray_cast.h"
  39.  
  40. double RayCast::CalcFormFactor( Vertex3 *pvertex, Instance
  41.     *penv )
  42. {
  43.   int i;            // Loop index
  44.   double ff;        // Vertex-source form factor
  45.   double ray_len;   // Ray length
  46.   Vector3 nv;       // Vertex normal
  47.   Vector3 n_ray;    // Normalized ray direction
  48.   Vector3 r_ray;    // Reverse normalized ray direction
  49.   Vector3 view;     // Source patch view vector
  50.  
  51.   start = Vector3(pvertex->GetPosn());
  52.   nv = pvertex->GetNormal();
  53.   view = start - src_center;
  54.  
  55.   // Determine whether source patch is backface
  56.   if (Dot(src_norm, view) < MIN_VALUE)
  57.     return 0.0;
  58.  
  59.   ff = 0.0;
  60.   for (i = 0; i < RC_NumRays; i++)
  61.   {
  62.     // Select random point on source patch
  63.     Select(&end);
  64.  
  65.     // Generate ray to shoot from vertex to source
  66.     ray_dir = end - start;
  67.  
  68.     // Check for source point behind vertex
  69.     if (Dot(nv, ray_dir) < MIN_VALUE)
  70.       continue;
  71.  
  72.     // Test for ray-element intersection
  73.     if (CheckOcclusion(penv) == FALSE)
  74.     {
  75.       // Calculate ray length
  76.       ray_len = ray_dir.Length();
  77.  
  78.       // Calculate normalized ray direction
  79.       n_ray = ray_dir;
  80.       n_ray.Norm();
  81.  
  82.       // Determine reverse normalized ray direction
  83.       r_ray = -n_ray;
  84.  
  85.       // Update form factor estimation
  86.       ff += Dot(n_ray, nv) * Dot(r_ray, src_norm) / ((PI *
  87.           ray_len * ray_len) + ray_area);
  88.     }
  89.   }
  90.  
  91.   // Multiply by ray-source patch intersection area
  92.   ff *= ray_area;
  93.  
  94.   return ff;
  95. }
  96.  
  97. // Initialize parameters for source patch
  98. void RayCast::Init( Patch3 *ppatch )
  99. {
  100.   double a1, a2;        // Triangle areas
  101.   Vector3 temp;         // Temporary 3-D vector
  102.   Vector3 e0, e1, e2;   // Edge vectors
  103.  
  104.   psrc = ppatch;
  105.   pcache = NULL;
  106.   src_area = psrc->GetArea();
  107.   src_norm = psrc->GetNormal();
  108.   src_center = Vector3(psrc->GetCenter());
  109.   ray_area = src_area / RC_NumRays;
  110.  
  111.   // Get patch vertex vectors
  112.   v0 = Vector3(ppatch->GetVertexPtr(0)->GetPosn());
  113.   v1 = Vector3(ppatch->GetVertexPtr(1)->GetPosn());
  114.   v2 = Vector3(ppatch->GetVertexPtr(2)->GetPosn());
  115.   v3 = Vector3(ppatch->GetVertexPtr(3)->GetPosn());
  116.  
  117.   // Calculate patch edge vectors
  118.   e0 = Vector3(v1 - v0);
  119.   e1 = Vector3(v2 - v0);
  120.  
  121.   // Calculate first triangle area
  122.   temp = Cross(e0, e1);
  123.   a1 = temp.Length() / 2.0;
  124.  
  125.   if (ppatch->IsQuad() == TRUE)
  126.   {
  127.     // Calculate patch edge vector
  128.     e2 = Vector3(v3 - v0);
  129.  
  130.     // Calculate second triangle area
  131.     temp = Cross(e1, e2);
  132.     a2 = temp.Length() / 2.0;
  133.   }
  134.   else
  135.     a2 = 0.0;
  136.  
  137.   // Calculate fractional area of first triangle
  138.   selector = a1 / (a1 + a2);
  139. }
  140.  
  141. // Select random point within source patch area
  142. void RayCast::Select( Vector3 *ppoint )
  143. {
  144.   double s, t;      // Random point parameters
  145.  
  146.   // Get random point parameters
  147.   s = GetNormRand();
  148.   t = GetNormRand();
  149.  
  150.   // Ensure random point is inside triangle
  151.   if (s + t > 1.0)
  152.   {
  153.     s = 1.0 - s;
  154.     t = 1.0 - t;
  155.   }
  156.  
  157.   // Calculate random point co-ordinates
  158.   if (GetNormRand() <= selector)
  159.   {
  160.     // Locate point in first triangle
  161.     *ppoint = (1.0 - s - t) * v0 + s * v1 + t * v2;
  162.   }
  163.   else
  164.   {
  165.     // Locate point in second triangle
  166.     *ppoint = (1.0 - s - t) * v0 + s * v2 + t * v3;
  167.   }
  168. }
  169.  
  170. // Check for ray occlusion
  171. BOOL RayCast::CheckOcclusion( Instance *pinst )
  172. {
  173.   Patch3 *ppatch;       // Patch pointer
  174.   Surface3 *psurf;      // Surface pointer
  175.  
  176.   // Test cached patch for ray-patch intersection
  177.   if (TestPatch(pcache) == TRUE)
  178.     return TRUE;
  179.  
  180.   // Walk the instance list
  181.   while (pinst != NULL)
  182.   {
  183.     // Walk the surface list
  184.     psurf = pinst->GetSurfPtr();
  185.     while (psurf != NULL)
  186.     {
  187.       // Walk the patch list
  188.       ppatch = psurf->GetPatchPtr();
  189.       while (ppatch != NULL)
  190.       {
  191.         if (ppatch != psrc)     // Ignore source patch
  192.         {
  193.           // Test for ray-patch intersection
  194.           if (TestPatch(ppatch) == TRUE)
  195.           {
  196.             // Cache occluding patch
  197.             pcache = ppatch;
  198.  
  199.             return TRUE;
  200.           }
  201.         }
  202.         ppatch = ppatch->GetNext();
  203.       }
  204.       psurf = psurf->GetNext();
  205.     }
  206.     pinst = pinst->GetNext();
  207.   }
  208.  
  209.   return FALSE;
  210. }
  211.  
  212. // Check for ray-patch intersection (Badouel's Algorithm)
  213. BOOL RayCast::TestPatch( Patch3 *ppatch )
  214. {
  215.   BOOL i_flag;          // Intersection flag
  216.   int i;                // Loop index
  217.   int i0, i1, i2;       // Projection plane axis indices
  218.   double alpha;         // Scaling parameter
  219.   double beta;          // Scaling parameter
  220.   double dist;          // Patch plane distance
  221.   double d, t;          // Temporary variables
  222.   double isect[3];      // Ray-patch intersection
  223.   double n_mag[3];      // Patch normal axis magnitudes
  224.   double vert[4][3];    // Patch vertices
  225.   double s0, s1, s2;    // Projected vector co-ordinates
  226.   double t0, t1, t2;    // Projected vector co-ordinates
  227.   Point3 *pvp;          // Vertex position pointer
  228.   Vector3 normal;       // Patch normal
  229.   Vector3 temp;         // Temporary 3-D vector
  230.  
  231.   // Check for valid patch
  232.   if (ppatch == NULL)
  233.     return FALSE;
  234.  
  235.   // Get patch normal
  236.   normal = ppatch->GetNormal();
  237.  
  238.   // Calculate divisor
  239.   d = Dot(normal, ray_dir);
  240.  
  241.   // Determine whether ray is parallel to patch
  242.   if (fabs(d) < MIN_VALUE)
  243.     return FALSE;
  244.  
  245.   // Calculate patch plane distance
  246.   temp = Vector3(ppatch->GetVertexPtr(0)->GetPosn());
  247.   dist = Dot(normal, temp);
  248.       
  249.   // Calculate ray hit time parameter
  250.   t = (dist - Dot(normal, start)) / d;
  251.  
  252.   // Check whether patch plane is behind receiver vertex or
  253.   // source patch point
  254.   //
  255.   // NOTE: MIN_VALUE offsets are required to prevent
  256.   //       interpretation of adjoining surface vertices as
  257.   //       occlusions
  258.   if (t < MIN_VALUE || t > (1.0 - MIN_VALUE))
  259.     return FALSE;
  260.  
  261.   // Calculate ray-patch plane intersection
  262.   temp = start + (ray_dir * t);
  263.  
  264.   // Get intersection axes
  265.   isect[0] = temp.GetX();
  266.   isect[1] = temp.GetY();
  267.   isect[2] = temp.GetZ();
  268.  
  269.   // Get patch normal axis magnitudes
  270.   n_mag[0] = fabs(normal.GetX());
  271.   n_mag[1] = fabs(normal.GetY());
  272.   n_mag[2] = fabs(normal.GetZ());
  273.  
  274.   // Get patch vertex axes
  275.   for (i = 0; i < ppatch->GetNumVert(); i++)
  276.   {
  277.     pvp = ppatch->GetVertexPtr(i)->GetPosnPtr();
  278.     vert[i][0] = pvp->GetX();
  279.     vert[i][1] = pvp->GetY();
  280.     vert[i][2] = pvp->GetZ();
  281.   }
  282.  
  283.   // Find patch normal dominant axis
  284.   if ((n_mag[0] > n_mag[1]) && (n_mag[0] > n_mag[2]))
  285.   {
  286.     i0 = 0; i1 = 1; i2 = 2;     // X-axis dominant
  287.   }
  288.   else if ((n_mag[1] > n_mag[0]) && (n_mag[1] > n_mag[2]))
  289.   {
  290.     i0 = 1; i1 = 0; i2 = 2;     // Y-axis dominant
  291.   }
  292.   else
  293.   {
  294.     i0 = 2; i1 = 0; i2 = 1;     // Z-axis dominant
  295.   }
  296.  
  297.   // Calculate projected vertex #0 co-ordinates
  298.   s0 = isect[i1] - vert[0][i1];
  299.   t0 = isect[i2] - vert[0][i2];
  300.  
  301.   // Check for intersection (consider quadrilateral as two
  302.   // adjacent triangles
  303.   i = 2;
  304.   i_flag = FALSE;
  305.   do
  306.   {
  307.     // Calculate projected vertex co-ordinates
  308.     s1 = vert[i - 1][i1] - vert[0][i1];
  309.     t1 = vert[i - 1][i2] - vert[0][i2];
  310.  
  311.     s2 = vert[i][i1] - vert[0][i1];
  312.     t2 = vert[i][i2] - vert[0][i2];
  313.  
  314.     // Determine vector scaling parameters
  315.     if (fabs(s1) < MIN_VALUE)   // Is s1 == 0 ?
  316.     {
  317.       beta = s0 / s2;
  318.       if ((beta >= 0.0) && (beta <= 1.0))
  319.       {
  320.         alpha = (t0 - beta * t2) / t1;
  321.         i_flag = ((alpha >= 0.0) && ((alpha + beta) <=
  322.             1.0));
  323.       }
  324.     }
  325.     else
  326.     {
  327.       beta = (s1 * t0 - s0 * t1) / (s1 * t2 - s2 * t1);
  328.       if ((beta >= 0.0) && (beta <= 1.0))
  329.       {
  330.         alpha = (s0 - beta * s2) / s1;
  331.  
  332.         // Test for intersection
  333.         i_flag = ((alpha >= 0.0) && ((alpha + beta) <=
  334.             1.0));
  335.       }
  336.     }
  337.     i++;    // Advance to next triangle (if any)
  338.   }
  339.   while (i_flag == FALSE && i < ppatch->GetNumVert());
  340.  
  341.   return i_flag;
  342. }
  343.  
  344.